Link del repositorio de github aquí
#antes de salir
save.image(file="objetos_actuales.Rdata")
#cuando te metes
load("OUTPUT/objetos_actuales.Rdata")
library(tidyjson)
library(rjson)
library(tidyverse)
library(dplyr)
library(ggplot2)
library(lubridate)
library(readr)
library(plotly)
En nuestro seminario vamos a estudiar las altas hospitalarias por ictus, en los diferentes hospitales de Castilla y León durante cuatro años, y como estas, están relacionadas con distintos factores de riesgo.
Lo primero de todo, explicar que los factor de riesgo son características o circunstancias que atentan contra el equilibrio, contra la salud, y que causan enfermedades o muerte. [senado1999factores]
Durante el desarrollo de nuestro seminario vamos a trabajar con datos de dos factores de riesgo, los cuales son: la calidad del aire y valorar si los centros de desintoxicación los podemos considerar como factor de riesgo o como factor de prevención de sufrir un ictus.
La calidad del aire es otro factor de riesgo del ictus, quizás menos conocido, pero estudios confirman que “Niveles más elevados de NO el día de ingreso conlleva un aumento del riesgo de mortalidad” [suarez2024influencia]
Objetivos principales:
1.Estudiar la cantidad de altas hospitalarias por ictus que hay en cada provincia de Castilla y León y ver como varían dependiendo del sexo y de la edad.
2.Analizar si a peor calidad del aire hay más riesgo de sufrir un ictus.
3.Estudiar la relación entre centros de ayuda a la desintoxicación por provincia y el número de ictus.
Los datos han sido recabados de Datos abiertos del gobierno de españa que hasta el momento contiene cerca de 90000 conjuntos de datos, distribuidos en un amplio número de temáticas.También supera las 550000 distribuciones, Además de 313 iniciativas de datos abiertos.
Descripción:
Contamos con 1 set de datos en tipo .csv. Este formato, (Comma Separated Values) son archivos de texto que en pricipio deberían estar con los caracteres separados por comas, aunque en algunos de nuestros datos se usa punto y coma en lugar de la coma.
Este es nuestro archivo .csv.: -calidad_aire
También contamos con dos sets de datos .json.(JavaScript Object Notation) es un formato para almacenar e intercambiar datos.
Estos son nuestros archivos .json.: -altas_hospitalarias_ictus -centros_servicios_saisde Este último proviene del Instituto Nacional de Estadística (INE).
altas_hospitalarias_ictus <- fromJSON(file="INPUT/data/altas_hospitalarias_ictus.json")
#altas_hospitalarias_ictus
calidad_aire <- read_csv("INPUT/data/calidad_aire.csv")
#calidad_aire
centros_servicios_saisde <- fromJSON(file = "INPUT/data/centros_servicios_saisde.json")
#centros_servicios_saisde
source("INPUT/functions/ictus_por_provincia.R")
Para resolver el primer objetivo vamos a:
altas_hospitalarias_ictus %>%
spread_all() %>%
gather_object %>%
json_types %>%
count(name, type)
## # A tibble: 4 × 3
## name type n
## <chr> <fct> <int>
## 1 datasetid string 20090
## 2 fields object 20090
## 3 record_timestamp string 20090
## 4 recordid string 20090
altas_hospitalarias_ictus %>%
enter_object(fields) %>%
gather_object %>%
spread_all %>%
json_types%>%
count(name, type)
## # A tibble: 11 × 3
## name type n
## <chr> <fct> <int>
## 1 ambito_de_procedencia string 20040
## 2 area string 20090
## 3 dia_de_la_semana_en_la_fecha_del_alta string 20090
## 4 dia_de_la_semana_en_la_fecha_del_ingreso string 20090
## 5 edad number 20090
## 6 fecha_de_alta string 20090
## 7 fecha_de_ingreso string 20090
## 8 hospital string 20090
## 9 provincia string 20090
## 10 sexo string 20090
## 11 zona_basica_de_salud_del_paciente string 20038
Dejamos la tabla entera tabulada, para trabajar a partir de ahora con este objeto:
altas_hosp_ictus_tab <- altas_hospitalarias_ictus%>%
spread_all()%>%
enter_object("fields")%>%
spread_all()
altas_hosp_ictus_tab
## # A tbl_json: 20,090 x 27 tibble with a "JSON" attribute
## ..JSON document.id datasetid recordid record_timestamp fields.fecha_de_ingr…¹
## <chr> <int> <chr> <chr> <chr> <chr>
## 1 "{\"f… 1 altas-ho… d642d9b… 2024-09-09T01:0… 2023-12-30T05:40:00+0…
## 2 "{\"f… 2 altas-ho… 1d47d7f… 2024-09-09T01:0… 2023-12-30T03:16:00+0…
## 3 "{\"f… 3 altas-ho… 95a004c… 2024-09-09T01:0… 2023-12-30T00:07:00+0…
## 4 "{\"f… 4 altas-ho… aab2391… 2024-09-09T01:0… 2023-12-29T18:29:00+0…
## 5 "{\"f… 5 altas-ho… 337420e… 2024-09-09T01:0… 2023-12-29T00:43:00+0…
## 6 "{\"f… 6 altas-ho… 32120d2… 2024-09-09T01:0… 2023-12-28T21:00:00+0…
## 7 "{\"f… 7 altas-ho… 4b61391… 2024-09-09T01:0… 2023-12-28T18:53:00+0…
## 8 "{\"f… 8 altas-ho… 9bb0dbc… 2024-09-09T01:0… 2023-12-28T15:25:00+0…
## 9 "{\"f… 9 altas-ho… 02c92ba… 2024-09-09T01:0… 2023-12-28T12:59:00+0…
## 10 "{\"f… 10 altas-ho… 3f05cd4… 2024-09-09T01:0… 2023-12-28T02:29:00+0…
## # ℹ 20,080 more rows
## # ℹ abbreviated name: ¹fields.fecha_de_ingreso
## # ℹ 21 more variables: fields.ambito_de_procedencia <chr>, fields.area <chr>,
## # fields.zona_basica_de_salud_del_paciente <chr>, fields.provincia <chr>,
## # fields.edad <dbl>, fields.sexo <chr>, fields.fecha_de_alta <chr>,
## # fields.dia_de_la_semana_en_la_fecha_del_ingreso <chr>,
## # fields.dia_de_la_semana_en_la_fecha_del_alta <chr>, …
Obtener el número de casos por provincias para empezar a relacionarlo con los siguientes objetivos
Provincias <- altas_hosp_ictus_tab%>%
rename(Provincia=fields.provincia)%>%
select(Provincia)
Provincias$Provincia <- ifelse(Provincias$Provincia == "Ávila", "Avila", Provincias$Provincia)
#Provincias
OBtener numero de ictus
prov_e_ictus <- Provincias%>%
group_by(Provincia)%>%
summarise(numero_ictus=n())
ictus_por_provincia(Provincias$Provincia)
## Avila Burgos León Palencia Salamanca Segovia Soria
## 1178 2788 3513 1064 3743 1142 601
## Valladolid Zamora
## 4575 1486
df_ictus_por_provincia <- data.frame(provincia=prov_e_ictus$Provincia,
n_ictus=prov_e_ictus$numero_ictus)
df_ictus_por_provincia
## provincia n_ictus
## 1 Avila 1178
## 2 Burgos 2788
## 3 León 3513
## 4 Palencia 1064
## 5 Salamanca 3743
## 6 Segovia 1142
## 7 Soria 601
## 8 Valladolid 4575
## 9 Zamora 1486
Calculamos el número de personas que han sufrido un ictus para cada edad, para ello obtenemos el atributo edad de cada una de las personas.Además lo haremos de dos maneras:
edades_total <- altas_hosp_ictus_tab %>%
rename(Edad = edad) %>%
select(Edad)
edades_total
## # A tbl_json: 20,090 x 2 tibble with a "JSON" attribute
## ..JSON Edad
## <chr> <dbl>
## 1 "{\"fecha_de_ingr..." 90
## 2 "{\"fecha_de_ingr..." 83
## 3 "{\"fecha_de_ingr..." 72
## 4 "{\"fecha_de_ingr..." 78
## 5 "{\"fecha_de_ingr..." 76
## 6 "{\"fecha_de_ingr..." 76
## 7 "{\"fecha_de_ingr..." 81
## 8 "{\"fecha_de_ingr..." 74
## 9 "{\"fecha_de_ingr..." 78
## 10 "{\"fecha_de_ingr..." 53
## # ℹ 20,080 more rows
Importamos función que cuenta el total para cada una de las edades:
Se trata de una función programada más desde un punto de vista de programación, más adelante mostraremos como hacerlo desde un enfoque en el manejo de datos como se requiere en la asignatura y por ende, la gran mayoría de código quedará implementado de la segunda manera.
source("INPUT/functions/total_por_edad.R")
resultados_edad <- total_por_edad(edades_total)
print(resultados_edad)
## Edad Conteo
## 1 1 3
## 2 2 0
## 3 3 0
## 4 4 0
## 5 5 1
## 6 6 0
## 7 7 1
## 8 8 3
## 9 9 5
## 10 10 0
## 11 11 0
## 12 12 1
## 13 13 2
## 14 14 4
## 15 15 2
## 16 16 1
## 17 17 2
## 18 18 1
## 19 19 2
## 20 20 3
## 21 21 6
## 22 22 3
## 23 23 2
## 24 24 2
## 25 25 6
## 26 26 2
## 27 27 4
## 28 28 9
## 29 29 5
## 30 30 7
## 31 31 12
## 32 32 8
## 33 33 10
## 34 34 9
## 35 35 15
## 36 36 19
## 37 37 21
## 38 38 28
## 39 39 20
## 40 40 33
## 41 41 36
## 42 42 53
## 43 43 46
## 44 44 53
## 45 45 57
## 46 46 75
## 47 47 79
## 48 48 102
## 49 49 97
## 50 50 106
## 51 51 111
## 52 52 141
## 53 53 143
## 54 54 165
## 55 55 190
## 56 56 180
## 57 57 190
## 58 58 201
## 59 59 269
## 60 60 211
## 61 61 282
## 62 62 283
## 63 63 344
## 64 64 338
## 65 65 322
## 66 66 360
## 67 67 393
## 68 68 347
## 69 69 411
## 70 70 412
## 71 71 421
## 72 72 432
## 73 73 456
## 74 74 518
## 75 75 503
## 76 76 614
## 77 77 619
## 78 78 629
## 79 79 587
## 80 80 661
## 81 81 560
## 82 82 612
## 83 83 665
## 84 84 668
## 85 85 729
## 86 86 732
## 87 87 641
## 88 88 711
## 89 89 590
## 90 90 514
## 91 91 456
## 92 92 411
## 93 93 288
## 94 94 248
## 95 95 181
## 96 96 130
## 97 97 85
## 98 98 67
## 99 99 45
## 100 100 22
## 101 101 16
## 102 102 9
## 103 103 8
## 104 104 6
ictus_por_edad <- edades_total%>%
group_by(Edad)%>%
summarise(num_ictus=n())
ictus_por_edad
## # A tibble: 99 × 2
## Edad num_ictus
## <dbl> <int>
## 1 0 7
## 2 1 3
## 3 5 1
## 4 7 1
## 5 8 3
## 6 9 5
## 7 12 1
## 8 13 2
## 9 14 4
## 10 15 2
## # ℹ 89 more rows
Realizamos una regresión polinómica:
regresion_edad_ictus <- lm(num_ictus~poly(Edad,3),data=ictus_por_edad)
summary(regresion_edad_ictus)
##
## Call:
## lm(formula = num_ictus ~ poly(Edad, 3), data = ictus_por_edad)
##
## Residuals:
## Min 1Q Median 3Q Max
## -227.10 -81.20 -21.57 79.80 279.25
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 202.93 10.82 18.759 < 2e-16 ***
## poly(Edad, 3)1 1449.03 107.64 13.462 < 2e-16 ***
## poly(Edad, 3)2 -527.86 107.64 -4.904 3.87e-06 ***
## poly(Edad, 3)3 -1368.52 107.64 -12.714 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 107.6 on 95 degrees of freedom
## Multiple R-squared: 0.7943, Adjusted R-squared: 0.7879
## F-statistic: 122.3 on 3 and 95 DF, p-value: < 2.2e-16
Hemos decidido realizar una regresión polinómica ya que una lineal no se ajusta exactamente lo que queremos y nos daba un R^2 muy bajo, en cambio la polinómica nos da un R^2 de 0.7994 lo cual nos indica que no es un mal modelo. Esto sucede porque en la regresión lineal debe crecer linealmente hasta el infinito pero a ciertas edades es biologicamente imposible que aumente ya que los humanos tenemos una vida limitada. En los resultados del modelor polinomial vemos que la edad es un factor estadísticamente siginificativo.
Gráfico que explica el modelo
ggplot(data=ictus_por_edad, aes(x=Edad,y=num_ictus))+
geom_point()+
geom_smooth()
Se ve claramente como a edades bajas el número deictus es muy reducido pero a partir de la barrera de los 45-50 empieza a subir drásticamente hasta los 75-88 que empieza a disminuir ya que a estas edades comienza a disminuir la tasa de supervivencia ya que las personas fallecen por otras patologías o por causas naturales.
Conteo por sexo:
sexototal <- altas_hosp_ictus_tab %>%
rename(Sexo = sexo) %>%
select(Sexo)
total_sexo <- sexototal$Sexo
Llamamos a la función que nos cuenta el número de hombres y de mujeres que han sufrido un ictus:
source("INPUT/functions/total_por_sexo.R")
resultados_sexo <- total_por_sexo(total_sexo)
print(resultados_sexo)
## Sexo Conteo
## 1 Hombre 11124
## 2 Mujer 8966
ictus_por_sexo<- sexototal%>%
group_by(Sexo)%>%
summarise(num_ictus=n())
ictus_por_sexo
## # A tibble: 2 × 2
## Sexo num_ictus
## <chr> <int>
## 1 Hombre 11124
## 2 Mujer 8966
Un total de 11124 sufrieron un ictus y el número de mujeres fue de 8966.
Mostramos de manera visual con un gráfico de barras sencillo la cantidad de hombres y mujeres que han sufrido un ictus:
ictus_por_sexo%>%
ggplot(mapping=aes(x=Sexo, y=num_ictus))+
geom_bar(stat="identity", fill="purple")+
labs(x="Sexo",y="Numero de ictus")
Realizamos un ANOVA: Para ello necesitamos conocer la media de la edad y su desviación estandar.
media_desviacion <- altas_hosp_ictus_tab %>%
group_by(sexo) %>%
summarise(
media = mean(edad, na.rm = TRUE),
desviacion_estandar = sd(edad, na.rm = TRUE)
)
media_desviacion
## # A tibble: 2 × 3
## sexo media desviacion_estandar
## <chr> <dbl> <dbl>
## 1 Hombre 73.3 13.0
## 2 Mujer 78.4 13.2
anova_resultado <- aov(edades_total$Edad ~ sexototal$Sexo, data = altas_hosp_ictus_tab)
summary(anova_resultado)
## Df Sum Sq Mean Sq F value Pr(>F)
## sexototal$Sexo 1 131314 131314 765.3 <2e-16 ***
## Residuals 20088 3446613 172
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Análisis estadístico: Tomamos como Test Hipótesis: H0: Medias iguales –> p-valor>0.05 H1: Medias distintas –>p-valor<0.05 Tras realizar el anova(analisis de la varianza) con la edad como variable numérica y el sexo como variable categórica vemos que el p-valor es muy bajo rechazamos la H0, lo que nos indica que el sexo es un factor estadisticamente significativo para los ictus.
Calculamos la media de los factores más destacados del aire, agrupados por fecha y provincia.
media_factores_aire <-calidad_aire%>%
rename(provincia=Provincia)%>%
mutate(Fechas = my(Fecha))%>%
mutate(mes= month(Fechas))%>%
mutate(anyo= year(Fechas))%>%
group_by(Fechas, provincia,mes,anyo)%>%
summarise(
NO=mean(`NO (ug/m3)`, na.rm = TRUE),
O3=mean(`O3 (ug/m3)`, na.rm = TRUE),
PM25=mean(`PM25 (ug/m3)`, na.rm = TRUE)
)%>%
arrange(Fechas)
media_factores_aire
## # A tibble: 459 × 7
## # Groups: Fechas, provincia, mes [459]
## Fechas provincia mes anyo NO O3 PM25
## <date> <chr> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 2019-10-01 Avila 10 2019 2 52 NaN
## 2 2019-10-01 Burgos 10 2019 4.67 49.8 4
## 3 2019-10-01 León 10 2019 2.87 44.2 NaN
## 4 2019-10-01 Palencia 10 2019 5.57 41.4 3.5
## 5 2019-10-01 Salamanca 10 2019 3 61 5
## 6 2019-10-01 Segovia 10 2019 2 57 NaN
## 7 2019-10-01 Soria 10 2019 8.5 44 NaN
## 8 2019-10-01 Valladolid 10 2019 11.5 38.9 11
## 9 2019-10-01 Zamora 10 2019 7 47 NaN
## 10 2019-11-01 Avila 11 2019 2 60 NaN
## # ℹ 449 more rows
Necesitamos cambiar el formato de la fecha, para así poder trabajar y agrupar datos que contienen fechas, ya que el paquete lubridate permite hacer lo que necesitamos con las fechas, en formato fecha.
#para poder hacer el estudio de manera correcta de la calidad del aire
altas_hosp_ictus_fecha_bien <- altas_hosp_ictus_tab %>%
mutate(Fecha = ymd_hms(fields.fecha_de_ingreso))%>%
mutate(mes= month(fields.fecha_de_ingreso))%>%
mutate(anyo = year(fields.fecha_de_ingreso))
altas_hosp_ictus_fecha_bien$provincia <- ifelse(altas_hosp_ictus_fecha_bien$provincia == "Ávila", "Avila",altas_hosp_ictus_fecha_bien$provincia)
altas_hosp_ictus_fecha_bien
## # A tbl_json: 20,090 x 30 tibble with a "JSON" attribute
## ..JSON document.id datasetid recordid record_timestamp fields.fecha_de_ingr…¹
## <chr> <int> <chr> <chr> <chr> <chr>
## 1 "{\"f… 1 altas-ho… d642d9b… 2024-09-09T01:0… 2023-12-30T05:40:00+0…
## 2 "{\"f… 2 altas-ho… 1d47d7f… 2024-09-09T01:0… 2023-12-30T03:16:00+0…
## 3 "{\"f… 3 altas-ho… 95a004c… 2024-09-09T01:0… 2023-12-30T00:07:00+0…
## 4 "{\"f… 4 altas-ho… aab2391… 2024-09-09T01:0… 2023-12-29T18:29:00+0…
## 5 "{\"f… 5 altas-ho… 337420e… 2024-09-09T01:0… 2023-12-29T00:43:00+0…
## 6 "{\"f… 6 altas-ho… 32120d2… 2024-09-09T01:0… 2023-12-28T21:00:00+0…
## 7 "{\"f… 7 altas-ho… 4b61391… 2024-09-09T01:0… 2023-12-28T18:53:00+0…
## 8 "{\"f… 8 altas-ho… 9bb0dbc… 2024-09-09T01:0… 2023-12-28T15:25:00+0…
## 9 "{\"f… 9 altas-ho… 02c92ba… 2024-09-09T01:0… 2023-12-28T12:59:00+0…
## 10 "{\"f… 10 altas-ho… 3f05cd4… 2024-09-09T01:0… 2023-12-28T02:29:00+0…
## # ℹ 20,080 more rows
## # ℹ abbreviated name: ¹fields.fecha_de_ingreso
## # ℹ 24 more variables: fields.ambito_de_procedencia <chr>, fields.area <chr>,
## # fields.zona_basica_de_salud_del_paciente <chr>, fields.provincia <chr>,
## # fields.edad <dbl>, fields.sexo <chr>, fields.fecha_de_alta <chr>,
## # fields.dia_de_la_semana_en_la_fecha_del_ingreso <chr>,
## # fields.dia_de_la_semana_en_la_fecha_del_alta <chr>, …
Juntamos ambas tablas con sus atributos comunes que nos son de interes (mes, año y provincia)
aire_con_ictus <- full_join(altas_hosp_ictus_fecha_bien, media_factores_aire, by = c("mes","anyo", "provincia"))
aire_con_ictus_resumen <- aire_con_ictus%>%
rename(ambito=fields.ambito_de_procedencia,dia_semana=fields.dia_de_la_semana_en_la_fecha_del_ingreso)%>%
select(Fechas,mes,anyo,edad,sexo,dia_semana,provincia,NO,O3,PM25,ambito)
aire_con_ictus_resumen
## # A tibble: 20,102 × 11
## Fechas mes anyo edad sexo dia_semana provincia NO O3 PM25
## <date> <dbl> <dbl> <dbl> <chr> <chr> <chr> <dbl> <dbl> <dbl>
## 1 2023-12-01 12 2023 90 Hombre Sábado Avila 4 38 NaN
## 2 2023-12-01 12 2023 83 Hombre Sábado Valladolid 12.5 31 10
## 3 2023-12-01 12 2023 72 Mujer Viernes León 6.88 34.5 9.5
## 4 2023-12-01 12 2023 78 Hombre Viernes Salamanca 6.67 49.5 6
## 5 2023-12-01 12 2023 76 Mujer Jueves Palencia 5.2 35.6 6
## 6 2023-12-01 12 2023 76 Hombre Jueves León 6.88 34.5 9.5
## 7 2023-12-01 12 2023 81 Hombre Jueves Zamora 9 35 NaN
## 8 2023-12-01 12 2023 74 Hombre Jueves Avila 4 38 NaN
## 9 2023-12-01 12 2023 78 Hombre Jueves Salamanca 6.67 49.5 6
## 10 2023-12-01 12 2023 53 Hombre Jueves Segovia 6 42 NaN
## # ℹ 20,092 more rows
## # ℹ 1 more variable: ambito <chr>
Calculamos la media de cada una de las sustancias (NO, O3 y PM25) agrupadas por provincia, año y mes. Además, contamos el número de ictus que hay para cada uno de ellos.
calcula_media_sustancias_por_provincia_y_meses <-aire_con_ictus_resumen %>%
group_by(provincia, anyo, mes) %>%
summarize(num_ictus=n(),
avg_NO = mean(NO, na.rm = TRUE),
avg_O3=mean(O3, na.rm = TRUE),
avg_PM25=mean(PM25, na.rm = TRUE))
calcula_media_sustancias_por_provincia_y_meses
## # A tibble: 459 × 7
## # Groups: provincia, anyo [45]
## provincia anyo mes num_ictus avg_NO avg_O3 avg_PM25
## <chr> <dbl> <dbl> <int> <dbl> <dbl> <dbl>
## 1 Avila 2019 10 1 2 52 NaN
## 2 Avila 2019 11 1 2 60 NaN
## 3 Avila 2019 12 6 6 43 NaN
## 4 Avila 2020 1 24 5 46 NaN
## 5 Avila 2020 2 17 5 42 NaN
## 6 Avila 2020 3 29 3 66 NaN
## 7 Avila 2020 4 19 2 65 NaN
## 8 Avila 2020 5 20 2 69 NaN
## 9 Avila 2020 6 37 3 72 NaN
## 10 Avila 2020 7 24 2 89 NaN
## # ℹ 449 more rows
Estudios confirman que “Niveles más elevados de NO el día de ingreso conlleva un aumento del riesgo de mortalidad debido a ictus” El gráfico que se muestra a continuación nos sirve para ver sí cuando la concentración de NO es mayor hay un mayor número de ictus.
Grafico de dispersión del NO:
calcula_media_sustancias_por_provincia_y_meses%>%
ggplot(aes(x=avg_NO, y=num_ictus))+
geom_point(aes(colour = factor(provincia)))+
geom_smooth(method = "lm", aes(colour = factor(provincia)))
Podemos observar que las lineas para cada provincia según el articulo(ref) deberían salir con pendiente positiva pero no es así. Creemos que este resultado se debe a un pequeño matiz entre el artículo y nuestros datos. Dicho matiz se trata de que en el artículo nos indica que concentraciones altas de NO aumentaban el riesgo de mortalidad mientras que en nuestros datos los pacientes que tenemos han recibido el alta. Así que como no tenemos número de fallecidos no lo podemos contrastar realmente con el artículo.
El ozono es uno de los principales contaminantes del aire. Tal y como afirma el profesor Shaowei Wu, de la universidad Xi´an Jiaotong (China): “estudios previos han sugerido que la contaminación por ozono daña el corazón y los vasos sanguíneos”.
Por ello vamos a realizar un gráfico de dispersión en el que veamos como varía el número de ictus dependiendo de como sea la media de O3 en cada provincia.
Gráfico de dispersión del 03:
calcula_media_sustancias_por_provincia_y_meses%>%
ggplot(aes(x=avg_O3, y=num_ictus))+
geom_point(aes(colour = factor(provincia)))+
geom_smooth(method = "lm", aes(colour = factor(provincia)))
En este caso, sí que podemos observar una pendiente positiva a mayor concentración de ozono por lo que podemos concluir que si es un factor que aumenta el numero de ictus.
Las hospitalizaciones aumentan en un día determinado entre un 0,5% y un 1% por cada aumento de 10 mg/m3 de PM2.5. PM2.5: Hace referencia al material particulado respirable presente en la atmósfera. Se puede encontrar más información sobre PM2.5 aqui
Por ello vamos a realizar un gráfico de dispersión para ver si cuando los niveles de PM2.5 son más elevados hay mayor número de ictus.
Grafico de dispersión del PM25:
calcula_media_sustancias_por_provincia_y_meses%>%
ggplot(aes(x=avg_PM25, y=num_ictus))+
geom_point(aes(colour = factor(provincia)))+
geom_smooth(method = "lm", aes(colour = factor(provincia)))
En este caso tenemos provincias que no tienen valores de PM2.5 (las que no tienen linea en el gráfico) pero en las que lo tenemos en todas se cumple que a mayor concentración mayor número de ictus.
Ahora vamos a realizar 3 gráficos, uno para cada sustancia
Gráfico que muestra la media de NO (en una escala de colores) por provincia y mes:
calcula_media_sustancias_por_provincia_y_meses %>%
distinct(provincia, .keep_all = TRUE) %>%
ggplot(aes(x = reorder(provincia,anyo), y = num_ictus)) +
geom_bar(stat = "identity", aes(fill = avg_NO)) +
scale_fill_gradient(low = "turquoise", high = "blue" , space = "Lab", na.value = "grey50", guide = "colourbar") +
labs(x = "", y = "Número ictus", fill = "Media de concentración de NO") +
theme_classic()
En esta gráfica aparecen representados los números de ictus por
provincia y la concentración de NO, la cual aparece pintada más oscura
cuando es mayor y más clara cuando es menor.
p <- calcula_media_sustancias_por_provincia_y_meses %>%
ggplot(aes(x = interaction(provincia, anyo), y = num_ictus, fill = avg_NO)) +
geom_bar(stat = "identity") +
scale_fill_gradient(low = "turquoise", high = "blue", space = "Lab", na.value = "grey50", guide = "colourbar") +
labs(x = "Provincia y Año", y = "Número de ictus", fill = "Media de concentración de NO") +
theme_classic() +
theme(axis.text.x = element_text(angle = 90, hjust = 1))
p
# Convierte el gráfico ggplot a plotly
ggplotly(p, tooltip = c("x", "y", "fill"))
Esta gráfica es una versión mejorada de la anterior, en la cual aparecen representados los ictus que ha habido en cada provincia y en cada año, cada barra tiene 12 partes diferenciadas (una por mes) en las de los años 2020 en adelante; en el año 2019 aparecen menos ya que solo contamos con datos de octubre en adelante. Las barras están coloreadas de distintas tonalidades, corriespondiendose las oscuras a mayores concentraciones de NO y las claras a menores concentraciones de NO. Además pasando el cursor por las barras, gracias a la tecnología del paquete plotly, podemos ver podemos ver con que año y provincia se corresponde esa barra, la media de concentración de NO en ese mes y también el número total de ictus que han ocurrido ese mes.
Gráfico que muestra la media de 03 (en una escala de colores) por provincia y mes:
calcula_media_sustancias_por_provincia_y_meses %>%
distinct(provincia, .keep_all = TRUE) %>%
ggplot(aes(x = reorder(provincia,anyo), y = num_ictus)) +
geom_bar(stat = "identity", aes(fill = avg_O3)) +
scale_fill_gradient(low = "turquoise", high = "blue", space = "Lab", na.value = "grey50", guide = "colourbar") +
labs(x = "", y = "Número ictus", fill = "Media de concentración de 03") +
theme_classic()
En esta gráfica aparecen representados los números de ictus por
provincia y la concentración de O3, la cual aparece pintada más oscura
cuando es mayor y más clara cuando es menor.
p1 <- calcula_media_sustancias_por_provincia_y_meses %>%
ggplot(aes(x = interaction(provincia, anyo), y = num_ictus, fill = avg_O3)) +
geom_bar(stat = "identity") +
scale_fill_gradient(low = "turquoise", high = "blue", space = "Lab", na.value = "grey50", guide = "colourbar") +
labs(x = "Provincia y Año", y = "Número de ictus", fill = "Media de concentración de O3") +
theme_classic() +
theme(axis.text.x = element_text(angle = 90, hjust = 1))
# Convierte el gráfico ggplot a plotly
ggplotly(p1, tooltip = c("x", "y", "fill"))
Esta gráfica es una versión mejorada de la anterior, en la cual aparecen representados los ictus que ha habido en cada provincia y en cada año, cada barra tiene 12 partes diferenciadas (una por mes) en las de los años 2020 en adelante; en el año 2019 aparecen menos ya que solo contamos con datos de octubre en adelante. Las barras están coloreadas de distintas tonalidades, corriespondiendose las oscuras a mayores concentraciones de O3y las claras a menores concentraciones de O3. Además pasando el cursor por las barras, gracias a la tecnología del paquete plotly, podemos ver podemos ver con que año y provincia se corresponde esa barra, la media de concentración de NO en ese mes y también el número total de ictus que han ocurrido ese mes.
Gráfico que muestra la media de PM2.5 (en una escala de colores) por provincia y mes:
calcula_media_sustancias_por_provincia_y_meses %>%
distinct(provincia, .keep_all = TRUE) %>%
ggplot(aes(x = reorder(provincia,anyo), y = num_ictus)) +
geom_bar(stat = "identity", aes(fill = avg_PM25)) +
scale_fill_gradient(low = "turquoise", high = "blue", space = "Lab", na.value = "grey50", guide = "colourbar") +
labs(x = "", y = "Número ictus", fill = "Media de concentración de PM2.5") +
theme_classic()
En esta gráfica aparecen representados los números de ictus por
provincia y la concentración de O3, la cual aparece pintada más oscura
cuando es mayor y más clara cuando es menor.
p2 <- calcula_media_sustancias_por_provincia_y_meses %>%
ggplot(aes(x = interaction(provincia, anyo), y = num_ictus, fill = avg_PM25)) +
geom_bar(stat = "identity") +
scale_fill_gradient(low = "turquoise", high = "blue", space = "Lab", na.value = "grey50", guide = "colourbar") +
labs(x = "Provincia y Año", y = "Número de ictus", fill = "Media de concentración de PM2,5") +
theme_classic() +
theme(axis.text.x = element_text(angle = 90, hjust = 1))
# Convierte el gráfico ggplot a plotly
ggplotly(p2, tooltip = c("x", "y", "fill"))
Esta gráfica es una versión mejorada de la anterior, en la cual aparecen representados los ictus que ha habido en cada provincia y en cada año, cada barra tiene 12 partes diferenciadas (una por mes) en las de los años 2020 en adelante; en el año 2019 aparecen menos ya que solo contamos con datos de octubre en adelante. Las barras están coloreadas de distintas tonalidades, corriespondiendose las oscuras a mayores concentraciones de O3y las claras a menores concentraciones de O3. Además pasando el cursor por las barras, gracias a la tecnología del paquete plotly, podemos ver podemos ver con que año y provincia se corresponde esa barra, la media de concentración de NO en ese mes y también el número total de ictus que han ocurrido ese mes. Cabe destacar que en este caso en concreto, hay muchos valores NaN, que son los que aparecen de color gris.
Como dato curioso, queremos ver si el día de la semana está relacionado con los ictus
numero_ictus_dia_semana <-aire_con_ictus_resumen %>%
group_by(dia_semana) %>%
summarize(num_ictus=n())
numero_ictus_dia_semana %>%
ggplot(mapping=aes(x = dia_semana, y = num_ictus)) +
geom_bar(stat = "identity", fill="orchid") +
labs(x = "Días de la semana", y="Número de ictus") +
theme_classic()
Por lo leído y descrito en uno de los artículos, niveles altos de NO el día de ingreso conlleva a mas riesgo de mortalidad, al ser nuestros datos de altas hospitalarias por ictus y no de muertes, lo que vamos a intentar ver es si los días que hay un mayor número de ingresos coinciden con niveles altos de NO.
#Importamos los datos de los centros de desintoxicación Visionamos los datos y vemos que con hacer un spread_all podemos obtener lo que necesitamos para nuestro trabajo.
spread_all(centros_servicios_saisde)
## # A tbl_json: 99 x 22 tibble with a "JSON" attribute
## ..JSON document.id datasetid recordid record_timestamp fields.posicion
## <chr> <int> <chr> <chr> <chr> <dbl>
## 1 "{\"dataseti… 1 centros-… a42164a… 2024-04-09T16:1… NA
## 2 "{\"dataseti… 2 centros-… 3f4b02d… 2024-04-09T16:1… NA
## 3 "{\"dataseti… 3 centros-… 4118619… 2024-04-09T16:1… NA
## 4 "{\"dataseti… 4 centros-… 7b8fef6… 2024-04-09T16:1… NA
## 5 "{\"dataseti… 5 centros-… aa67d45… 2024-04-09T16:1… NA
## 6 "{\"dataseti… 6 centros-… 536cdf6… 2024-04-09T16:1… NA
## 7 "{\"dataseti… 7 centros-… ba153c0… 2024-04-09T16:1… NA
## 8 "{\"dataseti… 8 centros-… 557ce5f… 2024-04-09T16:1… NA
## 9 "{\"dataseti… 9 centros-… 96ef383… 2024-04-09T16:1… NA
## 10 "{\"dataseti… 10 centros-… 1706963… 2024-04-09T16:1… NA
## # ℹ 89 more rows
## # ℹ 16 more variables: fields.e_mail <chr>, fields.provincia_abreviado <chr>,
## # fields.entidad_dependencia_institucional <chr>, fields.c_postal <chr>,
## # fields.tipo_de_recurso <chr>, fields.telefono <chr>,
## # fields.nivel_asistencial <chr>, fields.fax <chr>, fields.localidad <chr>,
## # fields.longitud <dbl>, fields.direccion <chr>, fields.provincia <chr>,
## # fields.latitud <dbl>, geometry.type <chr>, geometry.coordinates <dbl>, …
Centros <- centros_servicios_saisde %>%
spread_all() %>%
mutate(fields.provincia = ifelse(fields.provincia == "Ávila", "Avila", fields.provincia)) %>%
count(provincia = fields.provincia) %>%
as.data.frame()
Centros
## provincia n
## 1 Avila 6
## 2 Burgos 12
## 3 León 18
## 4 Palencia 13
## 5 Salamanca 16
## 6 Segovia 5
## 7 Soria 4
## 8 Valladolid 17
## 9 Zamora 8
ordenadas <- arrange(Centros,Centros$provincia)
ordenadas
## provincia n
## 1 Avila 6
## 2 Burgos 12
## 3 León 18
## 4 Palencia 13
## 5 Salamanca 16
## 6 Segovia 5
## 7 Soria 4
## 8 Valladolid 17
## 9 Zamora 8
#Importamos función que imprime los habitantes de castilla y león por provincias
source("INPUT/functions/poblacion_cyl.R")
poblacion_cyl
## provincia Poblacion_total Poblacion_extranjera
## 1 Avila 159764 12868
## 2 Burgos 357370 33150
## 3 León 448573 23398
## 4 Palencia 157787 8581
## 5 Salamanca 327089 18606
## 6 Segovia 155332 20015
## 7 Soria 89528 10363
## 8 Valladolid 521333 33352
## 9 Zamora 166927 7601
## 10 Castilla y León 2383703 167934
## Porc_pob_extranj_sobre_pob_total Num_munic_pob_extranj_mas_3_por_ciento
## 1 8.05 106
## 2 9.27 220
## 3 5.21 102
## 4 5.43 88
## 5 5.68 130
## 6 12.88 167
## 7 11.57 83
## 8 6.39 150
## 9 4.55 110
## 10 7.04 1156
str(poblacion_cyl)
## 'data.frame': 10 obs. of 5 variables:
## $ provincia : chr "Avila" "Burgos" "León" "Palencia" ...
## $ Poblacion_total : num 159764 357370 448573 157787 327089 ...
## $ Poblacion_extranjera : num 12868 33150 23398 8581 18606 ...
## $ Porc_pob_extranj_sobre_pob_total : num 8.05 9.27 5.21 5.43 5.68 ...
## $ Num_munic_pob_extranj_mas_3_por_ciento: num 106 220 102 88 130 ...
Para tener una referencia entre los centros de desintoxicación por provincia con el número de habitantes de esa provincia calcularemos cuál es el número de habitantes por centro; con el obetivo de ver realmente si hay más ictus en aquellas provincias donde el número de habitantes por centro es menor y el consumo de drogas es un factor significativo en los ictus o por el contrario no lo es.
pob_centros <- full_join(poblacion_cyl,ordenadas)
tabla <- pob_centros%>%
select(provincia,Poblacion_total,n)
source("INPUT/functions/divide.R")
hab_centro <- divide(tabla$Poblacion_total,tabla$n)
#faltaria indicar que con que provincia se corresponden esos datos.
tabla$hab_centro <- hab_centro
resultado <- tabla %>%
select(provincia, hab_centro)
print(resultado)
## provincia hab_centro
## 1 Avila 26627.33
## 2 Burgos 29780.83
## 3 León 24920.72
## 4 Palencia 12137.46
## 5 Salamanca 20443.06
## 6 Segovia 31066.40
## 7 Soria 22382.00
## 8 Valladolid 30666.65
## 9 Zamora 20865.88
## 10 Castilla y León NA
Lo lógico es esperar que en las provincias que cuentan con un menor número de habitantes, haya menor número de ictus y menor número de centros de desintoxicación.
numero_ictus <- ictus_por_provincia(Provincias)
df <- as.data.frame(numero_ictus)
pob_prov_y_total <- poblacion_cyl%>%
select(provincia, Poblacion_total)
pob_prov_y_total
## provincia Poblacion_total
## 1 Avila 159764
## 2 Burgos 357370
## 3 León 448573
## 4 Palencia 157787
## 5 Salamanca 327089
## 6 Segovia 155332
## 7 Soria 89528
## 8 Valladolid 521333
## 9 Zamora 166927
## 10 Castilla y León 2383703
pob_cyl_sin_ult_fila <- pob_prov_y_total [1:9,]
pob_cyl_sin_ult_fila
## provincia Poblacion_total
## 1 Avila 159764
## 2 Burgos 357370
## 3 León 448573
## 4 Palencia 157787
## 5 Salamanca 327089
## 6 Segovia 155332
## 7 Soria 89528
## 8 Valladolid 521333
## 9 Zamora 166927
arrange(.data=ordenadas,(ordenadas$provincia))
## provincia n
## 1 Avila 6
## 2 Burgos 12
## 3 León 18
## 4 Palencia 13
## 5 Salamanca 16
## 6 Segovia 5
## 7 Soria 4
## 8 Valladolid 17
## 9 Zamora 8
datos <- data.frame(pob=pob_cyl_sin_ult_fila,
centros= arrange(.data=ordenadas,(ordenadas$provincia)),
num_ictus=df$numero_ictus)
str(datos)
## 'data.frame': 9 obs. of 5 variables:
## $ pob.provincia : chr "Avila" "Burgos" "León" "Palencia" ...
## $ pob.Poblacion_total: num 159764 357370 448573 157787 327089 ...
## $ centros.provincia : chr "Avila" "Burgos" "León" "Palencia" ...
## $ centros.n : int 6 12 18 13 16 5 4 17 8
## $ num_ictus : num 1178 2788 3513 1064 3743 ...
datos
## pob.provincia pob.Poblacion_total centros.provincia centros.n num_ictus
## 1 Avila 159764 Avila 6 1178
## 2 Burgos 357370 Burgos 12 2788
## 3 León 448573 León 18 3513
## 4 Palencia 157787 Palencia 13 1064
## 5 Salamanca 327089 Salamanca 16 3743
## 6 Segovia 155332 Segovia 5 1142
## 7 Soria 89528 Soria 4 601
## 8 Valladolid 521333 Valladolid 17 4575
## 9 Zamora 166927 Zamora 8 1486
Gráfica en la que aparece representado el número de centros por provincia:
datos%>%
ggplot(mapping = aes(x = reorder(pob.provincia, num_ictus), y = centros.n)) +
geom_bar(stat = "identity",fill="blue") +
labs(x = "Provincias", y = "Número de centros") +
scale_y_continuous(expand = expansion(mult = c(0, .1))) +
theme_classic() +
theme(axis.text.x = element_text(angle = 60, vjust = 1, hjust = 1))
Gráfica de barras en la que se representa el número de ictus en cada provincia, coloreando de un color más fuerte las barras con los más centros de desintoxicación.
datos %>%
ggplot(mapping = aes(x = reorder(pob.provincia, num_ictus), y = num_ictus)) +
geom_bar(stat = "identity", aes(fill = centros.n)) +
scale_fill_gradient(low = "peachpuff", high = "red", space = "Lab", na.value = "grey50", guide = "colourbar") +
labs(x = "Provincias", y = "Número de ictus", fill = "Centros detox") +
scale_y_continuous(expand = expansion(mult = c(0, .1))) +
theme_classic() +
theme(axis.text.x = element_text(angle = 60, vjust = 1, hjust = 1))
Gráfica en la que representamos el número de habitantes por centro para cada provincia:
resultado
## provincia hab_centro
## 1 Avila 26627.33
## 2 Burgos 29780.83
## 3 León 24920.72
## 4 Palencia 12137.46
## 5 Salamanca 20443.06
## 6 Segovia 31066.40
## 7 Soria 22382.00
## 8 Valladolid 30666.65
## 9 Zamora 20865.88
## 10 Castilla y León NA
resultado_sin_ult_fila=resultado[1:9,]
str(resultado_sin_ult_fila)
## 'data.frame': 9 obs. of 2 variables:
## $ provincia : chr "Avila" "Burgos" "León" "Palencia" ...
## $ hab_centro: num 26627 29781 24921 12137 20443 ...
ggplot(mapping = aes(x = reorder(datos$pob.provincia, datos$num_ictus), y = resultado_sin_ult_fila$hab_centro)) +
geom_bar(stat = "identity", fill = "green") +
labs(x = "Provincias", y = "Número habitantes por centro") +
scale_y_continuous(expand = expansion(mult = c(0, .1))) +
theme_classic() +
theme(axis.text.x = element_text(angle = 60, vjust = 1, hjust = 1))
Obtenemos una tabla que nos permite visualizar de manera numérica toda la información estudiada anteriormente:
tabla_centros <- tibble(Provincia=datos$pob.provincia,
Numero_de_ictus=datos$num_ictus,
Poblacion_total_cyl=pob_cyl_sin_ult_fila$Poblacion_total,
Numero_centros=datos$centros.n,
Habitantes_por_centro=resultado_sin_ult_fila$hab_centro,
Porcentaje_ictus_por_hab=(Numero_de_ictus/Poblacion_total_cyl)*100)
tabla_centros
## # A tibble: 9 × 6
## Provincia Numero_de_ictus Poblacion_total_cyl Numero_centros
## <chr> <dbl> <dbl> <int>
## 1 Avila 1178 159764 6
## 2 Burgos 2788 357370 12
## 3 León 3513 448573 18
## 4 Palencia 1064 157787 13
## 5 Salamanca 3743 327089 16
## 6 Segovia 1142 155332 5
## 7 Soria 601 89528 4
## 8 Valladolid 4575 521333 17
## 9 Zamora 1486 166927 8
## # ℹ 2 more variables: Habitantes_por_centro <dbl>,
## # Porcentaje_ictus_por_hab <dbl>
Tras la realización de las tres gráficas anteriores para ver como están relacionados el número de ictus de cada provincia con el número de centros de desintoxicación que hay en las mismas, y el número de habitantes por centro de cada una de ellas, hemos llegado a las siguientes conclusiones:
Las poblaciones que cuentan con un menor número de habitantes (Ávila, Palencia, Soria y Zamora) son a su vez las que cuentan con un menor porcentaje de ictus; esto seguramente pueda deberse a que por regla general en las ciudades pequeñas se suelen llevar mejores habitos de vida.
Podemos destacar tres casos en especial:
En Palencia, cuentan con un número elevado de centros de desintoxicación para la población que tiene; lo cual supone que haya pocos habitantes por centro, pudiendo ver así como el porcentaje de ictus es el segundo más pequeño.Por lo tanto podría parecer que la presencia de más centros ayuda a reducir el número de ictus.
Segovia tiene un gran número de habitantes por centro, el más alto de todos, lo cual podría indicar que las personas que acuden van a recibir un “peor” tratamiento ya que el tratamiento se va a tener que repartir entre más pacietes y al ser peor el riesgo de ictus va a ser mayor, sin embargo en esta población no se cumple la hipotésis planteada en la conclusión de Palencia.
Salamanca es la provincia que cuenta con un mayor porcentaje de ictus por habitante, lo cual es destacable ya que no es la que más habitantes tiene y además es la segunda que menos habitantes por centro tiene.
Por lo tanto,visualmente y falta de un estudio estadístico preciso no podemos determinar si el número de centros es un factor estadísticamente significativo para el riesgo de padecer un ictus.
Queda demostrado que tanto la edad como el sexo influyen en el riesgo de sufrir un ictus, esto lo hemos comprobado realizando una regresión polinómica (edad) y un ANOVA (sexo), con sus correspondientes gráficas para hacerlo más visual.
En cuanto a las sustancias del aire, podemos decir que la concentración de NO y la de PM2,5 si que se relacionan con la probabilidad de sufrir un ictus; mientras que la concentración de O3 no podemos llegar a saber si verdaderamente esta relacionada ya que lo que está comprobado es que el hecho de que haya elevadas concentraciones de NO aumenta la mortalidad por ictus en ese día, y nuestros datos no son de muertes sino de altas hospitalarias.
En relación con la cantidad de centros, tampoco podemos determinar si el que haya más numero de centros reduce la probabilidad de sufrir un ictus.
.